Electronic correlation in nanoscale junctions: Comparison of the GW approximation 
to a numerically exact solution of the single-impurity Anderson model 
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The impact of electronic correlation in nanoscale junctions, e.g. formed by single molecules, 
is analyzed using the single-impurity Anderson model. Numerically exact Quantum Monte Carlo 
calculations are performed to map out the orbital filling, linear response conductance and spec- 
tral function as a function of the Coulomb interaction strength and the impurity level position. 
These numerical results form a benchmark against which approximate, but more broadly applica- 
ble, approaches to include electronic correlation in transport can be compared. As an example, the 
self consistent GW approximation has been implemented for the Anderson model and the results 
compared to this benchmark. For weak coupling or for level positions such that the impurity is 
either nearly empty or nearly full, the GW approximation is found to be accurate. However, for 
intermediate or strong coupling, the GW approximation does not properly represent the impact of 
spin or charge fluctuations. Neither the spectral function nor the linear response conductance are 
accurately given across the Coulomb blockade plateau and well into the mixed valence regimes. 
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I. INTRODUCTION 

Electrical circuits containing nanoscale junctions are 
central to nanoscience and condensed matter physics. 
The importance comes both from the conceptual issues 
they raise and the possibility of qualitatively smaller elec- 
tronic devices with new electrical properties 1]. Exam- 
ples include metallic atomic-scale point contacts that ex- 
hibit quantized conductance and organic molecules 
linking conducting leads. These latter may form non- 
resonant tunnel junctions 0,01 or single molecule devices 
whose conductance is controlled by a Kondo resonance 

siaa. 

The challenge to theory is twofold. First, the atomic 
scale specifics of chemical bonding and local structure can 
profoundly influence local potentials and energy align- 
ments for the electronic states that control conduction. 
Second, even if electron-electron interactions in the leads 
are well screened and may be effectively incorporated into 
the energy bands, interactions on the molecule are typi- 
cally not small, and may strongly affect the conductance 
and the spectrum. 

The important role played by the chemical and struc- 
tural details has led to a strong emphasis in the literature 
on self consistent theories, often utilizing approximate 
implementations of Density Functional Theory (DFT) 
In these theories the atomic scale potential and the 
local bonding structure are treated in detail while the 
conductance is calculated via a Landauer approach based 
on the electronic states derived from the self consistent 
Hamiltonian [§, However, these approaches treat 

excited states in a mean field manner and there remain 
significant questions concerning the role of electron cor- 
relations and whether the important energy levels are 
accurately represented by the mean field theories that 



are utiHzed 

The self consistent DFT approach has proved to be rel- 
atively accurate for metallic point contacts and certain 
molecular junctions with conductance near Go = 2e^//i, 
the quantum of conductance [H, [H, However, de- 
tailed comparison between experiment and self consis- 
tent calculations of conductance through single molecule 
junctions has generally shown a large discrepency, up to 
several orders of magnitude in the non-resonant tunneling 
regime P, [H, . Typically the measured conductance is 
smaller than the calculated conductance. Unfortunately, 
these comparisons are complicated by substantial vari- 
ability in the measured values for the same molecule and 
signficant uncertainty about the atomic scale structure 
of the junction near the single molecule link [sj. Cal- 
culations are typically performed for relatively idealized 
junction structures and the conductance can be sensi- 
tive to the local geometry for widely used thiol linkages 

The impact of uncertainties in the junction structure 
on the comparison (i.e. of not performing the calculations 
for the relevant bonding configuration) has recently been 
clarified, following the discovery that amine linked junc- 
tions produce single molecule junctions with reproducible 
conductance measurements, a result that was understood 
to derive from a selective bonding motif [l^ . A study of 
the impact of amine-gold link structures on junction con- 
ductance for benzenediamine gave strong support to the 
selectivity of the bonding and showed good agreement 
between theory and experiment for the distribution of 
conductance [20|. However, the magnitude of the cal- 
culated conductance exceeded the measured value by a 
factor of seven. This suggests that even after the junction 
structure is reliably accounted for, discrepancies remain 
and points to the importance of correlation effects be- 
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yond the the commonly used DFT based self consistent 
approach. 

These challenges have stimulated theoretical research 
along a number of lines, including implementing self in- 
teraction corrections [U, [2^ in order to obtain im- 
proved estimates of energy level alignment, correlated ba- 
sis function techniques (24 . [25j to improve the description 
of the electronic wave functions, fundamental analyses 
of the application of DFT to electronic transport to go 
beyond the Landauer approach while remaining within 
the DFT framework [2G, 27], and finally diagrammatic 
perturbation analyses of beyond-DFT correlation effects 

In this paper we focus on the last issue, namely the 
dynamical consequences of the on- molecule interactions. 
To discuss the essential physics, we restrict attention to 
a single resonance coupled to metallic leads and includ- 
ing the local Coulomb interaction, the single-impurity 
Anderson model [33j] . We do not consider many body ef- 
fects associated with interactions in the leads or between 
the molecule and the leads. The present work examines 
equilibrium properties and conductance in the linear re- 
sponse regime only. 

We present a numerically exact quantum Monte Carlo 
solution to this simple model as a benchmark against 
which other approximate approaches can be compared. 
As an example of such a comparison, we analyse the 
GW approximation for the electron self energy [U, . 
The GW approach [s^, HH] has been applied in full de- 
tail to succesfuUy predict the quasiparticle energies for a 
wide range of solids, surfaces, molecules and nanosystems 
[H, 113, S, m, Motivated by these successes, and 
noting that it is a conserving approximation [4l], \^ , sev- 
eral authors have begun to apply the GW approximation 
to calculate the electronic properties of single-molecule 
conductors [sol [sil . [3^ . However, it is far from clear 
under what circumstances the GW approximation will 
accurately treat the local correlations and the resulting 
transport phenomena. 

The rest of this paper is organized as follows: Section 
II presents the model, section III describes the methods 
(touching on the issue of the proper definition of the GW 
approximation for a local orbital), section IV presents re- 
sults for the level filling, conductance, and spectral func- 
tion, section V discuss the GW approximation and sec- 
tion VI is a conclusion. An Appendix presents some de- 
tails of the derivation of the GW equations we use. 



[331 , represented by the Hamiltonian 



H ^ ^Ekcl^Cka +'^£ddU'^ (1) 



Here dj. creates an electron of spin a on the localized level 
(energy Sd) and the U term describes the d ~ d interac- 
tion, cj.^ creates an electron of spin a in the lead state 
with energy Sk- Because we will be concerned only with 
equilibrium properties, a restriction to a single electronic 
continuum is possible: in a two-lead situation one combi- 
nation of lead states decouples from the problem and the 
state created by really refers to an electron in the ap- 
propriate "hybridizing" linear combination. Vk describes 
the hybridization between the level and the lead. 

The crucial quantity that describes the lead electrons 
is the hybridization function A(e) = tt^j, \ Vk\'^S{e — Sk)- 
In our work, we assume a semicircular density of states 
and a /c-independent V: 



A(e) = V- 



2^2 



lei < 2t. 



(2) 



We choose parameters so that A(e = 0) ^ t but our 
conclusions do not depend in any important way of this 
assumption. 

We shall be interested in correlations of the d-electrons, 
in particular the retarded Green's functionfi^, 

/"OO 

GdaH = -^ e»(-+»o^)*(K(i), 4(0)]), (3) 

Jo 

from which we obtain the spectral function (index "d" is 
dropped) : 



A„[uj) = -- ImG^(tj) 

TT 

The d-occupancy {n„) is given by 



dUilAcr{uo)f{Lo) 



(4) 



(5) 



Here the Fermi function /(w) = l/(exp(/3a;) -I- 1) and we 
have chosen the zero of energy such that the chemical 
potential jj, = 0. The linear response conductance cr is 
given by [4J]: 



II. MODEL 

We study the simplest possible model of a molecular 
junction: a single level which may hold 0, 1 or 2 elec- 
trons, has an interaction term which controls the energy 
of the two electron state, and is coupled to an electronic 
continuum which represents the leads and is taken to 
be noninteracting. This is the Anderson impurity model 



A^iio). (6) 



The non-interacting {U ~ 0) model can be solved 
exactly m SI; we obtain 
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LJ~Sd~^^{L0) 



(7) 
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with the lead self energy 



(8) 



For the semicircular density of states the lead self energy 
has a simple analytical form (although attention must be 
paid to the branch cut structure). We have, on the real 
and imaginary frequency axes respectively 



E^(w) = -iA{uj)Q{2t- 



+V 



2 uj - sgnwe(|w| - 2t)y/uj^ - 4^2 
2l2 



,(9) 



S^(iw„) 



^^2 »^« - zsgn(cj„)v/(u;„)2 + 4^2 ^^^^ 



For U ^ the model is no longer analytically solvable. 
The effect of the many-body interaction is expressed 
mathematically by the self energy {uj), defined by the 
relation 



G(^) 



1 



(11) 



We now qualitatively discuss the behaviour of the 
model. If we assume that the the hybridization is weak 
{V <C t) and the energy range we are considering is 
well inside the band {2t ^ \uj\), we can take A(e) — A 
{— V'^ /t in our case) and neglect the real part of so 
that = —iA. The important parameter is U/A. For 
U/ A — > the occupancy varies smoothly with Ed and the 
spectral function has a single, approximately Lorentzian, 
peak centered at Sd with half-width A: 
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TT (a; - erf)2 + A2 ' 



(12) 



For A ^ {U/A oo) we have an isolated ion de- 
coupled from the leads. There are four states: the empty 
state |0) with energy E = 0, the fully occupied state |ti) 
with energy E — 2ed + U and a magnetic doublet |t) or 
I J,) with energy E — Ed- The T — Q spectral function 
depends on the occupancy: If Ed > 0, the ground state 
is |0) and the spectral function consists of an addition 
peak at UJ — Ed- If < but Ed + U > the ground 
state is one of 1 1) or | J,) and the spectral function has 
a removal peak at uj = Ed < and an addition peak at 
UJ — Ed + U > 0. Finally, ii Ed + U < the ground state is 
Iti) and the spectral function has only a removal peak, 
centered at u) — Ed + U < 0. 

These elementary considerations suggest that (pro- 
vided the d-level occupancy is neither zero nor two) 
there exists a critical Uc at which the single-peaked spec- 
tral function characteristic of small U/A changes to the 
multi-peaked form found in the large U approximation. 
A reasonable estimate for the relevant t/-scale is pro- 
vided by the Hartree-Fock (HF) approximation [sl, lisj 
which, for the model used here yields Uc/ A = tt at 
occupancy n — I. In the Hartree-Fock approxima- 
tion the interaction term Un-\ni is approximated by 



^^("-t)"-! + ^"T("i) -U{n^){ni) implying = 
so that 



G.(c^) 



1 



Lu-{Ed + u{n^^))-j:yiij) 



(13) 



with (?i-ct) fixed from Eq.Q. 

The Hartree-Fock approximation incorrectly predicts 
that for U > Uc the ground state is spin polarized. 
Corrections to the Hartree-Fock approximation allow the 
spin to fluctuate, leading to the Kondo effect [46l. l47j . The 
ground state is non-magnetic, characterized by a Kondo 
energy scale given approximately bvjisj 



Tk « 0.2V2AC/ exp[7re<j(£d + U)/{2AU)]. 



(14) 



Equation (fH)) is vahd only if EdiEd + U) < 0. Quahta- 
tively, Eq. ([H]) shows that the Kondo temperature Tj. is 
minimal at the half filling point Ed — —U/2. 

The Kondo ground state is a Fermi liquid, for which 
the low frequency behavior of the many-body self energy 
is 

T^'^iuj) = U{n) + Eo + (1 - Z-^)uj + e{uj^, T^). (15) 

Here the U (n) is the Hartree shift in the d-level energy 
and Eq is any extra chemical potential shift arising from 
interactions beyond Hartree-Fock. An important conse- 
quence of Eq. (|15p is that at sufficiently low temperatures 



A{uj = 0) = - 



1 A{uj = 0) 



TT E 



*2 



+ A{UJ = 0)2 



(16) 



with £*^ = Ed + Re (E^(cj = 0) + E'^(cj = 0)) so that at 
density n = 1 {e'^ = 0), A{uj — 0) — and from Eq. (O, 
the conductance a 2e'^/h. 



III. METHODS 

In this section we describe both the GW method and 
the numerically exact Quantum Monte Carlo (QMC) 
method to which we compare it. 

A. GW 

In the GW approximation [s^, |3^ one defines a 
screened interaction W and approximates the electron 
self-energy as 



(iuJn) = -Ty^GcrjiuJn 



iv„,)W^{iv„^) (17) 



here written as a function of Matsubara frequency [4^ , T 
stands for temperature in unit of energy. In the extended 
solid state problem for which the GW approximation 
was originally introduced, W is taken to be the screened 
Coulomb interaction in the charge channel. In the im- 
purity model considered here it is essential to retain the 
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spin channel, which controls the low energy physics. Care 
must also be taken to respect the Pauli principle: the GW 
approximation corresponds to a partial resummation of 
the infinite set of diagrams which define the theory and 
one must ensure that this partial resummation includes 
all the diagrams necessary to respect antisymmetry. 
We rewrite the Hubbard interaction lining as a 2 x 2 

matrix in spin space V with components 



T4/3 = f/(l - 5a/3). 



(18) 



An alternative definition Vap — U is sometimes used in 
the literature. The two definitions are compared in Ap- 
pendix. 

The screened interaction is derived from the irreducible 
polarizability P through 



W = I-VP 



V. 



(19) 



In the GW approximation, no vertex corrections are in- 
cluded, so the polarizability is just the Random Phase 
Approximation bubble, 
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+ U/2 

FIG. 1: Comparison of GW, Hartree-Fock and QMC results 
for d-occupancy. Parameters: U = 1.05 (top), 2.1, 4.2, 8.4 
(bottom). V — 2.55, t = 10, A « 0.65. Left column: compar- 
ison of HF and GW (where (n^) and (n^) differ the two values 
are shown). Right column: comparison of QMC results to HF 
and GW results for average density per spin {n| -I- n|)/2. 



PaiiuJr, 



) = -T^GaiiuJm - iVn)Ga{iVn). 



Explicitly, the screened interaction is then 



U^P-^iQn) 



l-U^PA^a)P-A^n)' 



(20) 



(21) 



For later use we note that th e q uantity W defined in 
Eq. ([T9|l may be expressed [s^, as a correlation func- 
tion which for a paramagnetic ground state on the imag- 
inary time contour is 

W^Ar) - ^ [(rp(r)p(O)) + (ra(r)a(O))] (22) 



with p — n^+Tii — (n-^+iii) and a — 7i| — n^. Eq. (|22p will 
be used below in our analysis of the differences between 
the GW approximation and the exact results. 

Equations (fTTl) , (1^ , (PT|) define a self-consistent set of 
equations which are solved numerically by iteration. All 
quantities are calculated on a real frequency grid, with a 
frequency range of ±4t and a frequency spacing as small 
as t/W^ . Where convergence issues arise (intermedi- 
ate coupling and the non-magnetic phase) we use Pulay 
mixing [32, ISlj . 



B. QMC 

A numerically exact solution to the Anderson impu- 
rity model may be obtained using Quantum Monte Carlo 
techniques. For most of the results obtained here we used 
the Hirsch-Fye method [H, [s^; for some of the lowest 
temperature data we used the recently developed contin- 
uous time method [S^. [s^. The QMC calculations were 



mostly performed on a parallel computer cluster with 20 
dual core 2.2 GHz nodes; a typical point requires about 
10 hours of computer time per CPU. For Hirsch-Fye ei- 
ther 256 or 512 time slices were used and the lowest ac- 
cessible temperature was T = 0.025. Convergence was 
verified by comparing two different time slices or by com- 
parison to the continuous time method. For the continu- 
ous time method the perturbation orders were typically 
20 ~ 80, but at the lowest T orders up to ~ 140 were 
needed. We typically use 10"* time slices for the lowest 
accessible temperature T « 0.006. 

The central object in the calculation is the imaginary 
time Green's function, related to the spectral function 
via 



G(r) 



duj 



A(w)e-™ 

1 _|_ 



(23) 



d-electron density and spin correlation functions were 
also measured. Inversion of Eq. (|^5|) to obtain A{u!) from 
a computed G{t) is a numerically ill-posed problem. 
We used the maximum entropy method^56i]: while this 
method sometimes produces unphysical feature, no dif- 
ficulties were encountered in the results described here. 
Once A{uj) is determined, ReG'(a;) is obtained from the 
Kramers-Kronig relation and then the sclf-cncrgy from 
Eq. nil). 

It should be noted that the QMC method is formu- 
lated at T > 0. The computational expense increases 
rapidly as T — s- 0, limiting the temperatures which can 
be reached. 
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IV. RESULTS 
A. d-occupancy 

Fig. [1] shows the d-occupancy as a function of ed + U/2 
(the U/2 shift puts the particle- hole symmetric point at 
zero) . The left panel shows the level occupancy obtained 
from the HF and GW approximations. The unphysical 
magnetic solutions occurring at U > Uc arc visible as a 
difference between (n|) and (n|). (The QMC calcula- 
tions yield (rif) = (n^) at all {U,T) studied). For the 
parameters studied, U^^ ~ 2. The [/^^ depends more 
strongly on temperature than does making it dif- 
ficult to determine with precision. We find U^"^ k, 5 
at T = 0. The larger Uc for the GW approximation 
arises from the self-consistency (using G rather than t/o 
to compute P) ; the renormalization of the Green function 
suppresses the instability [s^, Hsl . 

The magnetic phase transition is an artifact of the HF 
and GW approximation schemes, but one may expect 
that the spin-averaged quantities are reasonably reliably 
represented. In the right column we plot the spin aver- 
aged (n). For strong interaction one can see clearly the 
Coulomb blockade plateau. All three methods yield a 
roughly correct shape for the occupancy vs level energy 
curve, in particular giving approximately correct widths 
for the Coulomb blockade plateau. At very weak inter- 
action strength, [U = 1.05) all three methods agree in 
detail. As the interaction is increased, differences ap- 
pear between the approximate and exact results. The 
differences are most pronounced near the edges of the 
Coulomb blockade plateau, in the mixed valence regime 
where charge fluctuations are significant. 



B. Spectral function 

The computed spectral function A{lij) is shown in 
Fig. H Focus first on the QMC resuhs over the wide 
frequency range (left panels, solid lines). At C/ = 1.05 w 
/2 the spectral function is very close to the non- 
interacting value. Moving to the middle panel, we see 



that when U is increased from U = 2.1 



HF 



to 



U — 4.2, the Hubbard bands begin to form, however 
a central peak remains. At T = the height of the cen- 
tral peak should be I/ttA « 0.5. The reduced height 
A(oj = 0) « 0.4 is an effect of the non-zero temper- 
ature used in the simulations. For U — 4.2, Eq. ITU) 
implies Tk ~ 0.04 approximately equal to the studied 
T. For U — 8.4 the Hubbard band is well formed, and 
the central peak now clearly interpreted as a "Kondo 
resonance" remains. It is interesting that traces of the 
Kondo resonance are visible even though the tempera- 
ture T — 0.05 studied is much greater than the Kondo 
temperature Tk « 0.004 estimated from Eq. (fTil) as pre- 
viously noted by Meir et al [11] . 

The dotted hues in Fig. [5] show the results of the GW 
approximation. The left panels show that GW agrees rea- 
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FIG. 2: Electron spectral function calculated by QMC (solid 
line) and GW (dashed line) at the particle-hole symmetric 
point ed + U/2 = for parameters T = 0.05, V = 2.55, t = 10, 
A « 0.65 with U = 1.05, 2.1, 4.2, 8.4 (top to bottom panels), 
Left panel: wide frequency range; right panel: expanded view 
of low frequency range. For U = 8.4, GW approximation pro- 
duces a magnetic solution; we show the spin-averaged spectral 
function. Due to the fact that A{uj) — A{—u!) at half-filling 
point, only positive frequencies are shown in the right panel. 
The non-interacting {U = 0) result is also shown as a dashed 
line in the top panels. We believe the very weak dip visible 
for u! near in the QMC calculation for [/ = 2.1 is an artifact 
of the analytical continuation procedure. It is not important 
for our subsequent discussion. 



sonably well with the exact results at J7 < 2.1 « U^^ , 
but at U = 4.2 < J/fw ^Yie GW does not produce the 
Hubbard bands and underestimates the hight of the cen- 
tral peak. At U = 8.4 > U^^ the GW approximation 
by contrast produces the Hubbard bands but misses the 
central peak. The right column is an expanded view of 
the central peak. For U = 1.05 ~ U^^/2 the two meth- 
ods agree well with each other, essentially because the 
interaction corrections are weak. As the correlations are 
increased, differences appear. We see that even in the 
U « U^^ regime where the GW approximation is reason- 
ably accurate, the low frequency lineshape is incorrect, 
with being too low near a; = and too high in the 

wings of the central peak. The differences become more 
severe for higher U. 



C. Conductance 

Figure [3] compares the QMC and GW predictions for 
the linear response conductance at several U values. We 
see that the GW approximation systematically underesti- 
mates the conductance, with noticeable differences from 
the QMC values even for the smallest JJ-value, U = 2.0 « 
U^^ , where the GW and QMC spectral functions agree 
reasonably well. We also note that general Fermi-liquid 
arguments imply that as T ^ 0, cr = <jt=o — T^/Q^ 
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U=0, T=0 
QMCT=0.1 '-e- 
GWT=0.1 




FIG. 3: Linear response conductance as a function of T cal- 
culated at half-filling point ed + U/2 = for different U, with 
V = 2.55, t = 10, A ~ 0.65. The Kondo temperatures es- 
timated from Eq. (O are ^ 0.09, 0.04, 0.004 respectively 
for U = 2.1,4.2,8.4. 



with Q a temperature scale of order T^. However, none 
of the calculations reveal a clear behavior except for 
the QMC calculations at [/ = 2.1; we expect that this is 
because in all of the other cases the Kondo temperature 
is close to or below the temperatures studied. 

The three panels of Fig. [5] show the dependence of a 
on level position Sd at two different temperatures. At 
T = we expect an approximately Lorentzian resonance 
lincshape, broadened from the non-interacting value by 
the density-dependent level-shift encoded in the real part 
of the self-energy. As T is increased the conductance 
decreases; the decrease from the T = value is a con- 
sequence of many-body scattering; It is expected to be 
most pronounced at the particle-hole symmetric point 
£d + U/2 = 0. This may be seen mathematically from 
Eq. (|14p for the Kondo temperature. In physical terms, 
the conductance involves valence fluctuation from the 
state n=ltort = 0orn = 2;at the half filled point, 
these states are most widely separated in energy, so the 
fluctuations are most easily disrupted by temperature. 

The top panel shows results for U = 2.1 « , along 
with the U = curve for comparison. The increased 
width of the interacting curve relative to the noninter- 
acting one is evident as is the approximately Lorenztian 
form. The GW and QMC results agree in the wings of 
the curve, but disagree in the small ed + U/2 regime, with 
the GW approximation overestimating the suppression of 
conductance by thermal fluctuations. 

The middle panel shows results for [/ = 4.2 « 2J7f ^ 
at two different temperatures. The QMC curves dis- 
play the theoretically expected evolution with temper- 
ature and level position. In the wings of the lineshape 



U=0, T=0 
QMC T=0.025 
QMC 1=0.1 '-e- 

GW T=0.025 

GW 1=0.1 




ed+U/2 



FIG. 4: Linear response conductance as a function of level 
position ed + U/2 calculated for U — 2.1 (upper panel) , U — 
4.2 and U = 8.4 (bottom panel) with V = 2.55, t = 10, A ^ 
0.65 at temperatures indicated. The noninteracting U — 0, 
T = result is also shown for comparison. 
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(say for \ed + U/2\ > 2) the curves are temperature in- 
dependent (for the temperatures studied) and have an 
approximately Lorentzian decay. In the central region 
(led + U/2\ < 1.5) the n{ed) curves shown in Fig. [1] in- 
dicate the beginning of a Coulomb blockade plateau and 
we see correspondingly a strongly temperature depen- 
dent suppression of the conductance. For these parame- 
ters the Kondo temperature estimated from Eq. (jl4p is 
« 0.04 for Ed + U/2 = 0; we see that for our lower tem- 
perature T = 0.025 « 0.6Tfe the conductance approaches 
the noninteracting value, as expected. While the GW 
approximation shows T dependence, it is too small for 
£d + U/2 = while extending strong T dependence too 
far with respect to the level position. The shape of the 
ed-dependencc of the conductance for a given T is gener- 
ally wrong through the Coulomb blockade region. 

Finally, the lowest panel shows results for the strongest 
couphng, U — 8.4 « AU^^ . Reference to Fig. [1] shows 
that for this interaction strength, the Coulomb block- 
ade plateau is well formed. The theoretically estimated 
Kondo temperature at Ed + U/2 = is « 0.004 rather 
lower than the lowest temperature studied; correspond- 
ingly the QMC conductance in the Coulomb blockade 
regime is small and strongly temperature dependent. In 
this regime the GW approximation predicts a magnetic 
state with a gap at the Fermi energy and no Kondo reso- 
nance, so that the conductance at small ed + U/2 is qual- 
itatively incorrect. The GW approximation produces the 
correct scale of ed + U/2 at which conductance is restored 
(because it produces a Coulomb blockade plateau of the 
correct width) but gives an incorrect description of the 
details of the conductance as a function of level position 
until beyond the edge of the plateau. 



D. Self energy 

Figure [5] shows the real frequency behavior of the 
self energy. We observe that the low frequency part of 
Rel](w) is linear with a negative slope, which is consis- 
tent with Eq. (dni), and the statement that Z < 1. QMC 
and GW give essentially the same values of Z for all in- 
teraction strengths shown. The weak structure visible 
for the smallest U near cj = is a numerical artifact of 
the analytical continuation. At higher frequencies differ- 
ences between GW and QMC are evident. In particular, 
the high frequency tail of the GW curve disagrees with 
the exact analytical result (see Fig. [7] below) . Turning 
now to the imaginary part of the self energy we first note 
that the low frequency part of ImS(a;) is approximately 
quadratic, which is consistent with Eq. (jl5p . Again the 
weak structures visible very close to w = are believed 
to be artifacts of the analytical continuation procedure. 
The non-zero value at w = is a temperature effect. 
The QMC lmT.{uj) is highly peaked; however GW fails 
to produce these peaks. 
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FIG. 5: Real (left column) and imaginary (right column) 
parts of real axis self energy calculated using QMC (solid 
lines) and GW (dashed lines) at the particle-hole symmetric 
point ed + U/2 = for parameters V = 2.55, f = 10, A « 0.65, 
T = 0.05 and [7=1.05 (top), 2.1, 4.2(bottom). The Hartree 
shift [/ (n) is subtracted from ReE. Note diflferent vertical 
axes on left and right panels. 
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FIG. 6: Comparison of the real part of screened interaction 
M^(iw„) from QMC (Eq. ^) and GW (Eq. (gTI) at Si + 
U/2 = 0. lmW{iu„) = 0. Inset: the ratio Wgw/W^qmc- 
Parameters: T = 0.1, U = 4.2, Ed = ~U /2, V = 2.55, t = 10, 
A ^ 0.65. 



V. ANALYSIS OF THE GW APPROXIMATION 

In this section the input to the GW approximation is 
compared to the QMC results. The results provide an ex- 
plicit measure of the relative importance of the neglected 
diagrams for the irreducible polarizability and the self 
energy. 

Figure[S]compares the screened interaction W obtained 
from the GW equation Eq. ([2T|) to the W obtained from 
the QMC calculation (Eq. ^) for [/ = 4.2 « 2Uf^ 
at the half filling point ed + U/2 = 0. We see that the 
interaction strongly increases as ujn 0; this is a sig- 
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FIG. 7: ?7-self-energy on Matsubara axis (upper panel) and 
real axis (lower panels) at half-filling Ed + U/2 = computed 
using QMC, GW and a hybrid scheme in which the self en- 
ergy is computed according to Eq. (|17|l with the QMC G and 
W. The analytic asymptotic behavior T,{iLj„) — / (iiujn) 
is also shown. . At half filling ReS(iLj„) — U/2 and U/2 is 
also subtracted from ReE(a;) shown in the lower left panel. 
Parameters: T = 0.1, U = 4.2, V = 2.55, t = 10, A « 0.65. 

nature of the slow spin fluctuations which are important 
to the physics. The charge fluctuations are suppressed 
by the Coulomb blockade effect. Remarkably, the exact 
and GW results agree very well at the lowest frequencies 
(see inset). However, the W from the GW approximation 
is too large at intermediate to large frequencies; screen- 
ing is incorrectly estimated. The ratio VFgw/W^qmc is 
also shown. The ratio quantifies the effect of vertex cor- 
rections, which are neglected in the GW approximation. 
The ratio becomes a constant above a certain frequency. 

As shown in the Appendix, the GW approximation 
gives an incorrect treatment of the high frequency tail 
of the self energy, implying a violation of the Pauli prin- 
ciple. To fix this one may consider replacing the self- 
consistently determined G and W in Eq. ^Tl} by the ex- 
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FIG. 8: Comparison of the spectral functions at half-filling 
point. GqmcM^qmc curve is the analytical continuation result 
from E = iGqmcVKqmc- Parameters: T = 0.1, U — 4.2, 
V = 2.55, f = 10, A « 0.65. 



act Green function and polarizibility. We have used our 
QMC simulations to measure G and W (from Eq. ((22l) ) 
and used the results to compute E from Eq. pT]) . Rep- 
resentative results are shown in Fig. [71 The upper panel 
shows the imaginary axis self energy. The GW curve is 
seen to have an incorrect asymptotic behavior. Using 
the QMC G and W produces a curve with the correct 
high frequency limit but with an incorrect low frequency 
behavior. A similar effect is seen in the lower panels, 
which display the analytically continued curves. Both 
GW and "GqmcT^QMc" curves substantially underesti- 
mate the frequency dependence of S. Additional insight 
comes from the spectral functions shown in Fig [8l We 
see that the "GqmcM^qmc" curve fails to reproduce the 
Hubbard band structure, and gives an incorrect magni- 
tude at low frequency. 

This comparison explicitly shows that the diagrams ne- 
glected in the GW approximation for the self energy are 
significant in this regime. The use of the Wqmc auto- 
matically includes the vertex corrections for the screened 
interaction. Approximations beyond GW also include 
vertex corrections explicitly in the expression for the self 
energy [H, . Figures highlight the role they play 
for U > Uf^. 

VI. CONCLUSION 

In this paper we used a numerically exact Quantum 
Monte Carlo method to obtain results for the density, 
conductance and spectral function of the single-impurity 
Anderson model, a simple theoretical paradigm problem 
for molecular conductors. The Anderson model is char- 
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acterized by two dimensionless combinations of three pa- 
rameters: a level position Ed, a level width A and an 
interaction U. In the conventions adopted in this paper, 
the impurity level is half filled when ed + U/2 = 0. A rel- 
evant measure of the interaction strength is provided by 
the Hartree-Fock approximation, which predicts a mag- 
netic state for interaction strengths greater than a crit- 
ical value. We define U^^ to be the critical value for 
the half filled level. For U > U^^ the impurity spectral 
function is characterized by a three-peak structure with 
upper and lower Hubbard bands at uj « ±11/ 2 and a 
central "Kondo" peak which controls the linear response 
conductance. 

Our results are intended as benchmarks against which 
other, more approximate but more widely applicable 
methods may be compared. We compared our results 
to those obtained from the GW approximation, a self- 
consistent partial resummation of diagrammatic pertur- 
bation theory. The GW method is attractive because 
it can be combined with band theory to yield material- 
specific results, but its efficacy at treating strong cor- 
relations is unclear. Recent literature has argued that 
the GW method provides a reasonable description of the 
physics of the low T (T = limit) transport proper- 
ties of molecular conductors for a range of intermediate 

We showed that for weak to moderate coupling regime 
{0 < U < U^^) or for a nearly full or nearly empty 
d-level, the GW approximation provides a reasonable de- 
scription. As the interaction U approaches U^^ , some 
systematic deviations are observed in the shape of the 
spectral function near a; = and the dependence of the 
conductance on the level position near Sd + U/2 = 0. For 
the intermediate coupling regime 

and the strong coupling regime {U > U^"^), the GW 
appproximation gives an inaccurate representation of the 
spectrum and the linear response conductance across the 
Coulomb blockade plateau. The accuracy is also limited 
for the mixed valence regions at the boundaries. 

These findings are summarized in Fig. [9l which shows 
the regions where GW does and does not work well in the 
plane of interaction strength and level position. "Works 
well" is of course an imprecise definition; in construct- 
ing Fig. [9] we defined "works well" as "GW conductance 
within 15% of QMC conductance at T = 0.1" . The crite- 
rion is temperature dependent as can be seen from Fig. [3] 
and interaction dependent, as can be seen from slope of 
the boundary line, which is less than 1/2. We find that 
the GW approximation is reliable when the level is tuned 
so that the density is far enough outside the Coulomb 
blockade region. For U = 4.2 the GW approximation be- 
comes reasonable for densities at the edge of the Coulomb 
blockade plateau, but for U ~ 8.4 the density must be 
tuned well away from the plateau before GW becomes 
accurate. 

Qualitatively, in the parameter regime in which GW 
produces a central peak in the spectral function, it does 
not produce the Hubbard sidebands, while the Hubbard 
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FIG. 9; Phase plane of interaction strength U and bare 
level energy Ed + U/2 shows regions where GW approxima- 
tion works poorly (white, region I) and well (grey, region II). 
Our criterion for the latter is \a^^ - ^r'^^^l/cr'^^^ < 15%. 
Parameters: T = 0.1, V = 2.55, t = 10, A ^ 0.65. 



sidebands are produced only as a consequence of an un- 
physical magnetic ordering instability. Papers in the 
literature interpret the central peak found in the non- 
magnetic GW calculation as a Kondo resonance. We 
believe this interpretation is not correct. It is interesting 
to note, however, that some aspects of the many body 
physics (for example the low frequency "mass renormal- 
ization" dT,/doj or the spin correlation function) are cor- 
rectly given by GW. This has been seen for the self con- 
sistent second order self energy as well [st} . 

The self consistent GW approximation has the virtue 
of being a conserving approximation [4ll . . In physical 
systems where the local molecular levels remain nearly 
filled or nearly empty, or where the hybridization is large, 
our results show that the GW approximation will be rea- 
sonably accurate. This suggests that an approach based 
on the GW approximation may be very useful for molec- 
ular conductors in the non-resonant tunneling regime 
where large discrepencies exist between theory and exper- 
iment. However, when the local Coulomb interactions on 
the molecule are strong, the GW approximation does not 
accurately represent the impact of local spin and charge 
fluctuations. Neither the spectral distribution nor the 
linear response conductance are given properly. Applica- 
tion of the GW approximation to nanoscale junctions in 
the Kondo regime is not well justified. The analysis of 
the screened Coulomb interaction W and the evaluation 
of the GW approximation for the self energy operator 
with the the exact (QMC) G and W showed that vertex 
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FIG. 10: Comparison of the spin-dependent and spin- 
independent approaches to the GW approximation. Parame- 
ters used: T = 0, f/ = 4.2, V = 2.55, i = 10, A « 0.65. 



corrections are quite significant in these cases. Unfortu- 
nately, while there are systematic guidelines for including 
vertex corrections properly so as to maintain a conserv- 
ing approximation, the resulting theory is substantially 
more complex jSoj . 

Our conclusions are based on linear response. Other 
situations, in particular the out-of-equilibrium Coulomb 
Blockade regime remain to be studied. 
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APPENDIX 

This Appendix discusses technical details of the GW 
calculations. The local interaction depends on electron 
spin. One may consider two forms: 



Vo^p = U 



"spin — independent" , (A-1) 
"spin — dependent" . (A-2) 



We now show that the asymptotic high-frequency behav- 
ior of the GW self energy implies that the spin-dependent 



interaction Eq. (jA-2p is more appropriate than the spin- 
independent interaction Eq. (|A-ip . 

It is convenient to separate out the Hartree and 
Fock terms, writing E = S^^ -|- Y,^^ . For the spin- 
independent interaction the Hartree term for spin a is 
U{n^ + ni) while the Fock term is —U{na-); for the spin 
dependent interaction the Fock term vanishes and the 
Hartree term is U{ng). In either case we have 



■t^^GW,' \ 



de ImW^^(e) 
Imiy™(t = 0) 



(A-3) 



where TO stands for time-ordered. In the exact per- 
turbation theory analysis, the screened interaction W is 
related to the spin-spin correlation function through the 
polarizability H [s^] : 

w^°{t) ^ V5{t) + yn™(<)i/ (A-4) 

where: 

H™(t) = -^{T{M) - (n.(t))]K(0) - (n..(0))]}) 

(A-5) 

In the spin-dependent case, the (1 — 5a, p) term in the in- 
teraction implies that Wa^rj involves only the correlator 
for the opposite spin, so one finds the following asymp- 
totic behavior: 



{/2 



(A-6) 



On the other hand, for the spin-independent interac- 
tion, all spin indices are involved and one obtains 



(w oo) 



( {nan„> ) - (n^) {n„, ) ) . 



a, a' 



Thus we see that the spin-dependent interaction repro- 
duces approximately the analytically known asymptotic 
behavior of the self-energy S^(iLi;„) — ^[(ng-)(l— (rig-))] 
whereas the spin-independent interaction does not. The 
asymptotic behavior is only approximately reproduced 
because GW cannot account correctly for (n?) = (n^). 
We show in the main text that using the exact W yields 
the correct asymptotic behavior of S but still does not 
produce an accurate approximation at general lo. The 
spin-independent case provides a much worse approxi- 
mation, which would be wrong even if the exact corre- 
lations functions were used. This is an indication that 
the spin-dependent 2-particle interaction is to be pref- 
ered over the spin-independent one in the context of the 
GW approximation applied to the Anderson model. 

Fig. [10] shows a comparison of spin-dependent and 
spin-independent GW at zero temperature. Comparison 
to the lowest temperature QMC result shows that the 
lineshape calculated from the spin-dependent GW ap- 
proximation is closer to the QMC lineshape than is the re- 
sult of the spin-independent calculation. This is because 
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the spin-dependent approach is free of self-interaction fluctuations, whereas the spin-independent approach ac- 
effects and it accounts for some of spin-spin quantum counts only for density-density quantum fluctuations. 
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